rm(list = ls())

#setwd("~/Dropbox/Energiewende_Project/All_Replication_Files")
setwd("C:/Users/chris/Dropbox/Energiewende_Project/All_Replication_Files")

#Loading in Packages ----------------------------------------------------------------------------------

library(dplyr)
library(sandwich)
library(lmtest)
library(margins)
library(ggplot2)
library(scales)
library(gridExtra)
library(readr)
library(stargazer)
library(reshape2) 
library(cregg)
library(kableExtra)
library(knitr)
library(modelsummary)
options(knitr.table.format = "latex")

####################################################################################################################################
############################################Candidate Conjoint Experiment###########################################################
####################################################################################################################################

##########Data Prep:########
data_cand <- read_csv("data.cand.final.csv", col_types = cols(.default = "?"))
weights <- read_csv("weights_data.csv")

data_cand <- merge(data_cand, weights, by=c("tic"), all.x=T)


data_cand$migration_policy <- relevel(as.factor(data_cand$migration_policy), 
                                      "Refugee Camps Abroad")

data_cand$cand_gender <- relevel(as.factor(data_cand$cand_gender),
                                 "female")

data_cand$social_policy <- relevel(as.factor(data_cand$social_policy),
                                   "No Pension Reform")

data_cand$energy_policy <- relevel(as.factor(data_cand$energy_policy),
                                   "Continue at Current Speed")

data_cand$cand_party <- relevel(as.factor(data_cand$cand_party),
                                "CDU/CSU")
##########Helper Functions:###########

#function to create plot data:
create.plot.data <- function(mes){
  
  plot.data <- as.data.frame(mes[2:14,1:2])
  
  plot.data[14:18,] <- 0
  
  plot.data$factor <- rownames(plot.data)
  
  plot.data$factor[14:18] <- c("cand_genderfemale", "cand_party CDU/CSU", "energy_policyContinue at Current Speed", 'migration_policyRefugee Camps Abroad',"social_policyNo Pension Reform")  
  
  plot.data$factor <- as.factor(plot.data$factor)
  
  #introduce factor levels for different variables:
  
  plot.data$factor <- factor(plot.data$factor, levels =c(levels(plot.data$factor),
                                                         "Gender",
                                                         "Party",
                                                         "Energy Policy",
                                                         "Migration Policy",
                                                         "Social Policy"
  ))
  
  plot.data$factor <- factor(plot.data$factor,
                             levels = levels(plot.data$factor)[
                               c(21,9,10,11, 12, #energy policy
                                 22, 15,13,14, #migration policy
                                 23, 17, 16, 18, #social policy
                                 20, 3, 4, 5, 6, 7, 8, #party
                                 19, 1, 2 #gender
                               ) 
                             ])
  plot.data$factor <- factor(plot.data$factor, levels = rev(levels(plot.data$factor)))
  plot.data$factor
  
  plot.data$confintlower <- plot.data$Estimate - 1.96 * plot.data$`Std. Error`
  plot.data$confintupper <- plot.data$Estimate + 1.96 * plot.data$`Std. Error`
  return(plot.data)
}


#helper functions for running regressions

runreg_candidate <- function(dataset, subgroup){
  
  dats <- dataset[dataset[,subgroup] == 1,]
  
  r1 <- lm(cand_choice ~ energy_policy + social_policy + 
             migration_policy + cand_gender + cand_party,
           dats, weights=weights) 
  
  vc1 <- vcovCL(r1, cluster = dats$ResponseId)
  
  out <- list(coeftest(r1, vc1),
              length(r1$fitted.values))
  return(out)
}

runreg_candidate_cov<- function(dataset, subgroup){
  dats <- dataset[dataset[,subgroup] == 1,]
  
  xvars <- "energy_policy + social_policy + migration_policy + cand_gender + cand_party + region + age + gender + incomequintiles"
  yvar <- "cand_choice"
  fmla.reg = as.formula(paste0(yvar, '~', xvars))
  
  regout <- lm(fmla.reg,data = dats, weights=weights) 
  
  vcreg <- vcovCL(regout, cluster = dats$ResponseId)
  
  regrob <- coeftest(regout, vcreg)
  
  out <- list(regrob,
              length(regout$fitted.values))
  
  return(out)
}

runreg_candidate_int<- function(dataset, intervvariable, subgroup){
  dats <- dataset[dataset[,subgroup] == 1,]
  
  xvars <- paste0("energy_policy + social_policy + migration_policy + cand_gender + cand_party +", intervvariable," + energy_policy*",  intervvariable)
  yvar <- "cand_choice"
  fmla.reg = as.formula(paste0(yvar, '~', xvars))
  
  regout <- lm(fmla.reg,data = dats, weights=weights) 
  
  vcreg <- vcovCL(regout, cluster = dats$ResponseId)
  
  regrob <- coeftest(regout, vcreg)
  
  out <- list(regrob,
              length(regout$fitted.values))
  
  return(out)
}

runreg_candidate_cov_int<- function(dataset, intervvariable, subgroup){
  dats <- dataset[dataset[,subgroup] == 1,]
  
  xvars <- paste0("energy_policy + social_policy + migration_policy + cand_gender + cand_party +", intervvariable," + energy_policy*",  intervvariable,"+ region + age + gender + incomequintiles")
  yvar <- "cand_choice"
  fmla.reg = as.formula(paste0(yvar, '~', xvars))
  
  regout <- lm(fmla.reg,data = dats, weights=weights) 
  
  vcreg <- vcovCL(regout, cluster = dats$ResponseId)
  
  regrob <- coeftest(regout, vcreg)
  
  out <- list(regrob,
              length(regout$fitted.values))
  
  return(out)
}

runreg_candidate_robust <- function(dataset, intervvariable, subgroup){
  dats <- dataset[dataset[,subgroup] == 1,]
  
  xvars <- paste0("energy_policy + social_policy + migration_policy + cand_gender + cand_party +", intervvariable," + energy_policy*",  intervvariable," + climatehumans * energy_policy+ climate_policyn*energy_policy + left_right_4 *energy_policy + party_2021_2*energy_policy ")
  yvar <- "cand_choice"
  fmla.reg = as.formula(paste0(yvar, '~', xvars))
  
  regout <- lm(fmla.reg,data = dats, weights=weights) 
  
  vcreg <- vcovCL(regout, cluster = dats$ResponseId)
  
  regrob <- coeftest(regout, vcreg)
  
  out <- list(regrob,
              length(regout$fitted.values))
  
  return(out)
}


supportinteractions <- c("supportintervention", "highsupportintervention", 
                         "supportinterventionfrail", "highsupportinterventionfrail",
                         "supportinterventionecon", "highsupportinterventionecon",
                         "supportinterventionworkers", "highsupportinterventionworkers") 

#function to plot average marginal component effects 

margeffplot_cand <- function(res, n.sample , n.res, colors, titleadd = NULL, fileadd = NULL, legendtitle){
  
  plot.data <- lapply(res, create.plot.data)
  
  pl <- ggplot(plot.data[[1]])+
    geom_point(aes(x = Estimate, y = factor,  color = paste(colors[1], " n = ", n.sample[[1]], sep=""), shape = paste(colors[1], " n = ", n.sample[[1]], sep="")))+
    geom_segment(aes(x=confintlower, xend = confintupper,
                     y = factor, yend = factor, color = paste(colors[1], " n = ", n.sample[[1]], sep="")))+
    geom_vline(xintercept = 0, lty = "dotted")+
    theme_bw()+
    scale_y_discrete(NULL, drop=FALSE,
                     labels = rev(c("Energy Policy",
                                    "Continue at Current Speed",
                                    "Slow Down Energy Transition",
                                    "Subsidies for Companies",
                                    "Subsidies for Low Income Households",
                                    "Migration Policy",
                                    "Refugee Camps Abroad",
                                    "Integrate Refugees Better",
                                    "Keep Refugees Away",
                                    "Social Policy",
                                    "No Pension Reform",
                                    "Early Retirement OK",
                                    "Raise Pension Age",
                                    "Party",
                                    "CDU/CSU",
                                    "AfD",
                                    "Bündnis90/Die Grünen",
                                    "Die Linke",
                                    "FDP",
                                    "SPD",
                                    "Gender",
                                    "Female",
                                    "Male")))
  
  pl <- pl + theme(axis.text.y = element_text(
    face = rev(c('bold', 'plain', 'plain', 'plain','plain',
                 'bold', 'plain', 'plain', 'plain',
                 'bold', 'plain', 'plain', 'plain',
                 'bold', 'plain', 'plain', 'plain','plain', 'plain', 'plain',
                 'bold', 'plain', 'plain'))
  ))
  
  
  pl <- pl + ggtitle(paste("Treatment Effects Candidate Conjoint", titleadd, sep="")) + scale_color_discrete(name = "subgroups")+
    guides(color=guide_legend(title=legendtitle),
           shape = guide_legend(title=legendtitle))
  
  
  if(n.res > 1){
    pl.temp <- as.data.frame(plot.data[[2]])
    pl <- pl + 
      geom_point(data = pl.temp, aes(x = Estimate, y = factor, color = paste(colors[2], " n = ", n.sample[[2]], sep=""), shape = paste(colors[2], " n = ", n.sample[[2]], sep="")), position = position_nudge(y = -0.2))+
      geom_segment(data = pl.temp, aes( x=confintlower, xend = confintupper,
                                        y = factor, yend = factor, color = paste(colors[2], " n = ", n.sample[[2]], sep="")), position = position_nudge(y = -0.2))  
  }
  
  if(n.res > 2){
    
    pl.temp <- as.data.frame(plot.data[[3]])
    pl <- pl + 
      geom_point(data = pl.temp, aes(x = Estimate, y = factor, color = paste(colors[3], " n = ", n.sample[3], sep=""),  shape =paste(colors[3], " n = ", n.sample[3], sep="")), position = position_nudge(y = +0.2))+
      geom_segment(data = pl.temp, aes( x=confintlower, xend = confintupper,
                                        y = factor, yend = factor, color = paste(colors[3], " n = ", n.sample[3], sep="")), position = position_nudge(y = +0.2))  
  }
  
  if(n.res > 3){
    
    pl.temp <- as.data.frame(plot.data[[4]])
    pl <- pl + 
      geom_point(data = pl.temp, aes(x = Estimate, y = factor, color = paste(colors[4], " n = ", n.sample[4], sep=""),  shape =paste(colors[4], " n = ", n.sample[4], sep="")), position = position_nudge(y = -0.1))+
      geom_segment(data = pl.temp, aes( x=confintlower, xend = confintupper,
                                        y = factor, yend = factor, color = paste(colors[4], " n = ", n.sample[4], sep="")),position = position_nudge(y = -0.1))  
  }
  
  if(n.res > 4){
    
    pl.temp <- as.data.frame(plot.data[[5]])
    pl <- pl + 
      geom_point(data = pl.temp, aes(x = Estimate, y = factor, color = paste(colors[5], " n = ", n.sample[5], sep=""), shape =paste(colors[5], " n = ", n.sample[5], sep="")), position = position_nudge(y = +0.1))+
      geom_segment(data = pl.temp, aes( x=confintlower, xend = confintupper,
                                        y = factor, yend = factor, color = paste(colors[5], " n = ", n.sample[5], sep="")), position = position_nudge(y = +0.1))  
  }
  
  if(n.res > 5){
    
    pl.temp <- as.data.frame(plot.data[[6]])
    pl <- pl + 
      geom_point(data = pl.temp, aes(x = Estimate, y = factor, color = paste(colors[6], " n = ", n.sample[6], sep=""), shape =paste(colors[6], " n = ", n.sample[6], sep="")),position = position_nudge(y = -0.3))+
      geom_segment(data = pl.temp, aes( x=confintlower, xend = confintupper,
                                        y = factor, yend = factor, color = paste(colors[6], " n = ", n.sample[6], sep="")),position = position_nudge(y = -0.3))  
  }
  
  pl <- pl +theme(text = element_text(size=16),
                  plot.title = element_text(size=18, hjust = 0.5),
                  axis.text.y = element_text(size=16),
                  legend.position = "bottom",
                  legend.title = element_text(size=16), #change legend title font size
                  legend.text = element_text(size=16),
                  axis.title= element_text(size=16))+
    xlab("Average Marginal Component Effects on Support for Candidate")+
    scale_color_grey(start = 0.1, end = 0.5)
  
  
  if(n.res == 1){
    
    pl <- pl + 
      scale_color_manual(values=c("black"))+theme(legend.position = "none")
  }
  
  
  ggsave(filename = paste("Graphs/candidate_conjoint.weighted.", fileadd, ".png", sep=""),
         plot = pl,
         width = 12, height = 10)
  
  if(n.res>3){
    ggsave(filename = paste("Graphs/candidate_conjoint.weighted.", fileadd, ".png", sep=""),
           plot = pl,
           width = 12, height = 16)
  }
  print(pl)
  print( paste("Graphs/candidate_conjoint.weighted.", fileadd, ".png", sep=""))
  return(pl)
}


##########Base-line Regression:############

baselinereg <- runreg_candidate(data_cand, "all")

plot.data <- list(baselinereg[[1]])

weighted1 <- margeffplot_cand(plot.data , nrow(data_cand)/10 , 1, "", titleadd = "\n Weighted Regressions", fileadd = NULL, "")

#########split sample plot#################

plot.data.i <- list()
n.subsets <- list()
list.regs <- list()
subsets <- c("highsupportintervention", "lowsupportintervention")

for(i in 1:length(subsets)){
  data <- runreg_candidate(data_cand, subsets[i])
  plot.data.i[[i]] <- data[[1]]
  n.subsets[[i]] <- data[[2]]/10
}


weighted2<- margeffplot_cand(plot.data.i , n.subsets ,2, c("High", "Low"), ", by Belief in State Intervention, \n Weighted Regressions", "govtint", "Belief in State Intervention")+theme(axis.text.y = element_blank())

weighted_results <- ggpubr::ggarrange(weighted1, weighted2, nrow=1, legend="bottom", widths = c(0.7, 0.5))

ggsave(filename = "Graphs/candidate_conjoint.weighted_graph.png",
       plot = weighted_results,
       width = 20, height = 10)

##########Creating Regression Tables:###############
header.tr <- "DV: Support for Candidate"

#table 1: baseline and interaction with Belief in State Intervention
reg.1 <- runreg_candidate(data_cand, "all")

reg.2 <- runreg_candidate_cov(data_cand, "all")

reg.3 <- runreg_candidate_int(data_cand, "highsupportintervention", "all")

reg.4 <- runreg_candidate_cov_int(data_cand, "highsupportintervention", "all")

reg.5 <- runreg_candidate_int(data_cand, "supportintervention", "all")

reg.6 <- runreg_candidate_cov_int(data_cand, "supportintervention", "all")

row1 <- t(c("Demographic Controls",  "Yes",  "Yes", "Yes")) %>% as.data.frame()

table.1 <- modelsummary(list( reg.2[[1]], reg.4[[1]], reg.6[[1]]), stars=T,
                        gof_omit = "AIC|BIC",
                        coef_rename = c("energy_policySlow Down Energy Transition"= "Energy Transition Slow Down",
                                        "energy_policySubsidies for Competitiveness"= "Energy Transition Company Subsidies",
                                        "energy_policySubsidies for Low Income Households"= 	"Energy Transition Low Income H. Subsidies",
                                        "social_policyEarly Retirement OK"= 	"Early Retirement OK",
                                        "social_policyRaise Pension Age"= "Raise Pension Age",
                                        "migration_policyIntegrate Refugees Better"= "Integrate Refugees Better",
                                        "migration_policyKeep Refugees Away" = "Keep Refugees Away",
                                        "cand_gendermale" = "Candidate Male",
                                        "cand_partyAfD" = "AfD",
                                        "cand_partyBündnis90/ Grünen" = "Bündnis90/ Grünen",
                                        "cand_partyDie Linke"= "Die Linke",
                                        "cand_partyFDP" = "FDP",
                                        "cand_partySPD" = "SPD" ,
                                        "highsupportintervention" = "Belief in St. Int. Bin",
                                        "supportintervention"= "Belief in St. Int. Cont." ),
                        coef_omit = "region|age|gender|income", add_rows = row1,output = "latex")





tableheader.tr <- c(c(" " = 1,  header.tr = 3 ))

# set vector names 
names(tableheader.tr) <- c(" ", header.tr)



regtable.t1 <- table.1 %>% kableExtra::add_header_above(header=tableheader.tr) 


save_kable(regtable.t1, "Tables/candidate_conjoint.weighted.tex")



#table 2, interactions with different subsets of government policy:

list.table2 <- list()

for(inter in supportinteractions[3:8]){
  list.table2[[inter]] <- runreg_candidate_cov_int(data_cand, inter, "all")
}

list_gs_tab1 <- list.table2[["highsupportinterventionecon"]]
list_gs_tab2 <- list.table2[["highsupportinterventionfrail"]]
list_gs_tab3 <- list.table2[["highsupportinterventionworkers"]]


table2 <- modelsummary(list(list_gs_tab1[[1]], list_gs_tab2[[1]], list_gs_tab3[[1]]), stars=T,
                       gof_omit = "AIC|BIC",
                       coef_rename = c("energy_policySlow Down Energy Transition"= "Energy Transition Slow Down",
                                       "energy_policySubsidies for Competitiveness"= "Energy Transition Company Subsidies",
                                       "energy_policySubsidies for Low Income Households"= 	"Energy Transition Low Income H. Subsidies",
                                       "social_policyEarly Retirement OK"= 	"Early Retirement OK",
                                       "social_policyRaise Pension Age"= "Raise Pension Age",
                                       "migration_policyIntegrate Refugees Better"= "Integrate Refugees Better",
                                       "migration_policyKeep Refugees Away" = "Keep Refugees Away",
                                       "cand_gendermale" = "Candidate Male",
                                       "cand_partyAfD" = "AfD",
                                       "cand_partyBündnis90/ Grünen" = "Bündnis90/ Grünen",
                                       "cand_partyDie Linke"= "Die Linke",
                                       "cand_partyFDP" = "FDP",
                                       "cand_partySPD" = "SPD" ,
                                       "highsupportinterventionecon" = "Belief in St. Int. Bin. Econ",
                                       "supportinterventionecon"= "Belief in St. Int. Cont. Econ",
                                       "highsupportinterventionfrail" = "Belief in St. Int. Bin. Frail",
                                       "supportinterventionfrail" = "Belief in St. Int. cont. Frail",
                                       "highsupportinterventionworkers" = "Belief in St. Int. Bin. Workers",
                                       "supportinterventionworkers" = "Belief in St. Int. cont. Workers"
                       ),
                       coef_omit = "region|age|gender|income", add_rows = row1
                       ,output = "latex"
)


header.t2 <- "DV: Support for Candidate"


tableheader.t2 <- c(c(" " = 1,  header.t2 = 6 ))

names(tableheader.t2) <- c(" ", header.t2)


header.t3.1 <- "Just Economy"
header.t3.2 <- "Just Sick and Elderly"
header.t3.3 <- "Just Workers"


tableheader.t3 <- c(c(" " = 1,  header.t3.1 = 1,
                      header.t3.2 = 1,
                      header.t3.3=1))
names(tableheader.t3) <- c("Belief State Intervention", header.t3.1, header.t3.2, header.t3.3)


# set vector names 



regtable.t2 <- table2%>% kableExtra::add_header_above(header=tableheader.t3)  %>% kableExtra::add_header_above(header=tableheader.t2) 

save_kable(regtable.t2, "Tables/candidate_conjoint_diff_gs.weighted.tex")


#figures by type of intervention

plot.data.1 <- list()
n.subsets <- list()
list.regs <- list()
subsets1 <- c("highsupportinterventionecon", "lowsupportinterventionecon")

for(i in 1:length(subsets1)){
  data <- runreg_candidate_cov(data_cand, subsets1[i])
  plot.data.1[[i]] <- data[[1]]
  n.subsets[[i]] <- data[[2]]/10
}

int.1 <- margeffplot_cand(plot.data.1 , n.subsets ,2, c("High", "Low"), ",\n Belief in State Intervention, \nJust Economy", "govtint1", "Belief in State Intervention")


plot.data.2 <- list()
n.subsets <- list()
list.regs <- list()
subsets2 <- c("highsupportinterventionfrail", "lowsupportinterventionfrail")

for(i in 1:length(subsets2)){
  data <- runreg_candidate(data_cand, subsets2[i])
  plot.data.2[[i]] <- data[[1]]
  n.subsets[[i]] <- data[[2]]/10
}

int.2 <- margeffplot_cand(plot.data.2 , n.subsets ,2, c("High", "Low"), ",\n Belief in State Intervention, \n Just Sick and Elderly", "govtint2", "Belief in State Intervention")+theme(axis.text.y = element_blank())


plot.data.3 <- list()
n.subsets <- list()
list.regs <- list()
subsets3 <- c("highsupportinterventionworkers", "lowsupportinterventionworkers")

for(i in 1:length(subsets3)){
  data <- runreg_candidate(data_cand, subsets3[i])
  plot.data.3[[i]] <- data[[1]]
  n.subsets[[i]] <- data[[2]]/10
}

int.3  <- margeffplot_cand(plot.data.3 , n.subsets ,2, c("High", "Low"), ",\n Belief in State Intervention, \n Just Workers", "govtint3", "Belief in State Intervention")+theme(axis.text.y = element_blank())

intcand.all <- ggpubr::ggarrange(int.1, int.2, int.3, nrow=1, common.legend = TRUE, legend="bottom", widths = c(0.8, 0.5,0.5))

ggsave(filename = "Graphs/candidate_conjoint.indintsplit.weighted.png",
       plot = intcand.all,
       width = 22, height = 10)





#Table 1 repeated for only people who passed the mock vignette
reg.1v <- runreg_candidate(data_cand, "candidateattention")

reg.2v <- runreg_candidate_cov(data_cand, "candidateattention")

reg.3v <- runreg_candidate_int(data_cand, "highsupportintervention", "candidateattention")

reg.4v <- runreg_candidate_cov_int(data_cand, "highsupportintervention", "candidateattention")

reg.5v <- runreg_candidate_int(data_cand, "supportintervention", "candidateattention")

reg.6v <- runreg_candidate_cov_int(data_cand, "supportintervention", "candidateattention")

table3 <- modelsummary(list(reg.2v[[1]],reg.4v[[1]],reg.6v[[1]]),stars=T,
                       gof_omit = "AIC|BIC",
                       coef_rename = c("energy_policySlow Down Energy Transition"= "Energy Transition Slow Down",
                                       "energy_policySubsidies for Competitiveness"= "Energy Transition Company Subsidies",
                                       "energy_policySubsidies for Low Income Households"= 	"Energy Transition Low Income H. Subsidies",
                                       "social_policyEarly Retirement OK"= 	"Early Retirement OK",
                                       "social_policyRaise Pension Age"= "Raise Pension Age",
                                       "migration_policyIntegrate Refugees Better"= "Integrate Refugees Better",
                                       "migration_policyKeep Refugees Away" = "Keep Refugees Away",
                                       "cand_gendermale" = "Candidate Male",
                                       "cand_partyAfD" = "AfD",
                                       "cand_partyBündnis90/ Grünen" = "Bündnis90/ Grünen",
                                       "cand_partyDie Linke"= "Die Linke",
                                       "cand_partyFDP" = "FDP",
                                       "cand_partySPD" = "SPD" ,
                                       "highsupportintervention" = "Belief in St. Int. Bin",
                                       "supportintervention"= "Belief in St. Int. Cont." ),
                       add_rows = row1,output = "latex", coef_omit = "region|age|gender|income"
)


regtable.t3 <- table3 %>% kableExtra::add_header_above(header=tableheader.tr) 


save_kable(regtable.t3, "Tables/candidate_conjoint_attention.weighted.tex")



table3f <- modelsummary(list(reg.1v[[1]], reg.2v[[1]], reg.3v[[1]],reg.4v[[1]],reg.5v[[1]],reg.6v[[1]]),stars=T,
                        gof_omit = "AIC|BIC",
                        coef_rename = c("energy_policySlow Down Energy Transition"= "Energy Transition Slow Down",
                                        "energy_policySubsidies for Competitiveness"= "Energy Transition Company Subsidies",
                                        "energy_policySubsidies for Low Income Households"= 	"Energy Transition Low Income H. Subsidies",
                                        "social_policyEarly Retirement OK"= 	"Early Retirement OK",
                                        "social_policyRaise Pension Age"= "Raise Pension Age",
                                        "migration_policyIntegrate Refugees Better"= "Integrate Refugees Better",
                                        "migration_policyKeep Refugees Away" = "Keep Refugees Away",
                                        "cand_gendermale" = "Candidate Male",
                                        "cand_partyAfD" = "AfD",
                                        "cand_partyBündnis90/ Grünen" = "Bündnis90/ Grünen",
                                        "cand_partyDie Linke"= "Die Linke",
                                        "cand_partyFDP" = "FDP",
                                        "cand_partySPD" = "SPD" ,
                                        "highsupportintervention" = "Belief in St. Int. Bin",
                                        "supportintervention"= "Belief in St. Int. Cont." ),
                        output = "latex"
)


regtable.t3f <- table3f %>% kableExtra::add_header_above(header=tableheader.tr) 


save_kable(regtable.t3f, "Tables/candidate_conjoint_attention_full.weighted.tex")


plot.data.a <- list()
n.subsets <- list()
list.regs <- list()
subsets <- c("highsupportintervention", "lowsupportintervention")
data_att <- data_cand %>% filter(candidateattention==1)


for(i in 1:length(subsets)){
  data <- runreg_candidate_cov(data_att, subsets[i])
  plot.data.a[[i]] <- data[[1]]
  n.subsets[[i]] <- data[[2]]/10
}

margeffplot_cand(plot.data.a , n.subsets ,2, c("High", "Low"), " by Belief in State Intervention, \n Passed Attention Check Only", "govtint_attention", "Belief in State Intervention")


###regressions with additional interaction terms


robustregs <- list()

for(inter in supportinteractions){
  robustregs[[inter]] <-  runreg_candidate_robust(data_cand, inter, "all")
}


table.rf1 <- modelsummary(list(robustregs[["highsupportintervention"]][[1]]), stars=T,
                          gof_omit = "AIC|BIC",
                          coef_rename = c("energy_policySlow Down Energy Transition"= "Energy Transition Slow Down",
                                          "energy_policySubsidies for Competitiveness"= "Energy Transition Company Subsidies",
                                          "energy_policySubsidies for Low Income Households"= 	"Energy Transition Low Income H. Subsidies",
                                          "social_policyEarly Retirement OK"= 	"Early Retirement OK",
                                          "social_policyRaise Pension Age"= "Raise Pension Age",
                                          "migration_policyIntegrate Refugees Better"= "Integrate Refugees Better",
                                          "migration_policyKeep Refugees Away" = "Keep Refugees Away",
                                          "cand_gendermale" = "Candidate Male",
                                          "cand_partyAfD" = "AfD",
                                          "cand_partyBündnis90/ Grünen" = "Bündnis90/ Grünen",
                                          "cand_partyDie Linke"= "Die Linke",
                                          "cand_partyFDP" = "FDP",
                                          "cand_partySPD" = "SPD" ,
                                          "highsupportintervention" = "Belief in St. Int. Bin",
                                          "supportintervention"= "Belief in St. Int. Cont." ,
                                          "climatehumans" = "Belief in MM CC",
                                          "Climate_policyn" = "Support Climate Pol.",
                                          "party_2021_2" = "2021 Party Vote"),
                          estimate  = "{estimate} ({std.error})",
                          statistic = NULL,
                          output="latex")



tableheader.tr1 <- c(c(" " = 1,  header.tr = 1 ))

# set vector names 
names(tableheader.tr1) <- c(" ", header.tr)


regtable.trf1 <- table.rf1 %>% kableExtra::add_header_above(header=tableheader.tr1) 

save_kable(regtable.trf1, "Tables/candidate_conjoint_robust_full1.weighted.tex")



table.rf2 <- modelsummary(list(robustregs[["supportintervention"]][[1]]), stars=T,
                          gof_omit = "AIC|BIC",
                          coef_rename = c("energy_policySlow Down Energy Transition"= "Energy Transition Slow Down",
                                          "energy_policySubsidies for Competitiveness"= "Energy Transition Company Subsidies",
                                          "energy_policySubsidies for Low Income Households"= 	"Energy Transition Low Income H. Subsidies",
                                          "social_policyEarly Retirement OK"= 	"Early Retirement OK",
                                          "social_policyRaise Pension Age"= "Raise Pension Age",
                                          "migration_policyIntegrate Refugees Better"= "Integrate Refugees Better",
                                          "migration_policyKeep Refugees Away" = "Keep Refugees Away",
                                          "cand_gendermale" = "Candidate Male",
                                          "cand_partyAfD" = "AfD",
                                          "cand_partyBündnis90/ Grünen" = "Bündnis90/ Grünen",
                                          "cand_partyDie Linke"= "Die Linke",
                                          "cand_partyFDP" = "FDP",
                                          "cand_partySPD" = "SPD" ,
                                          "highsupportintervention" = "Belief in St. Int. Bin",
                                          "supportintervention"= "Belief in St. Int. Cont." ,
                                          "climatehumans" = "Belief in MM CC",
                                          "Climate_policyn" = "Support Climate Pol.",
                                          "party_2021_2" = "2021 Party Vote"),
                          estimate  = "{estimate} ({std.error})",
                          statistic = NULL,
                          output="latex")




regtable.trf2 <- table.rf2 %>% kableExtra::add_header_above(header=tableheader.tr1) 

save_kable(regtable.trf2, "Tables/candidate_conjoint_robust_full2.weighted.tex")



##################################Marginal Means Plots####################
mms <- cj(data_cand,   cand_choice~energy_policy + social_policy + migration_policy + cand_gender + cand_party, id = ~ResponseId,
          estimate = "mm")

mms

#function to create plot data:
create.plot.data.mm <- function(mms){
  
  plot.data <- mms %>% select(estimate, lower, upper)
  
  plot.data$factor <- c("Continue at Current Speed", "Slow Down Energy Transition", 
                        "Subsidies for Companies", "Subsidies for Low Income Households",
                        "No Pension Reform", "Early Retirement OK", "Raise Pension Age",
                        "Refugee Camps Abroad", "Integrate Refugees Better",
                        "Keep Refugees Away", "Female", "Male", "CDU/CSU", 
                        "AfD", "Bündnis90/Die Grünen", "Die Linke", "FDP", "SPD")
  
  plot.data$factor <- as.factor(plot.data$factor)
  
  #introduce factor levels for different variables:
  
  plot.data$factor <- factor(plot.data$factor, levels =c(levels(plot.data$factor),
                                                         "Gender",
                                                         "Party",
                                                         "Energy Policy",
                                                         "Migration Policy",
                                                         "Social Policy"
  ))
  
  plot.data$factor <- factor(plot.data$factor,
                             levels = levels(plot.data$factor)[
                               c(21,4,15,17, 18, #energy policy
                                 22, 14,9,10, #migration policy
                                 23, 12, 6, 13, #social policy
                                 20, 3, 1, 2, 5, 7, 16, #party
                                 19, 8, 11 #gender
                               ) 
                             ])
  plot.data$factor <- factor(plot.data$factor, levels = rev(levels(plot.data$factor)))
  plot.data$factor
  
  return(plot.data)
}

plot.data.mm <- create.plot.data.mm(mms)

pl.mm <- ggplot(plot.data.mm)+
  geom_point(aes(x = estimate, y = factor))+
  geom_segment(aes(x=lower, xend = upper,
                   y = factor, yend = factor))+
  geom_vline(xintercept = 0.5, lty = "dotted")+
  theme_bw()+
  scale_y_discrete(NULL, drop=FALSE,
                   labels = rev(c("Energy Policy",
                                  "Continue at Current Speed",
                                  "Slow Down Energy Transition",
                                  "Subsidies for Companies",
                                  "Subsidies for Low Income Households",
                                  "Migration Policy",
                                  "Refugee Camps Abroad",
                                  "Integrate Refugees Better",
                                  "Keep Refugees Away",
                                  "Social Policy",
                                  "No Pension Reform",
                                  "Early Retirement OK",
                                  "Raise Pension Age",
                                  "Party",
                                  "CDU/CSU",
                                  "AfD",
                                  "Bündnis90/Die Grünen",
                                  "Die Linke",
                                  "FDP",
                                  "SPD",
                                  "Gender",
                                  "Female",
                                  "Male")))+ 
  ggtitle(paste("Marginal Means Candidate Conjoint,
All Participants"))+
  theme(axis.text.y = element_text(
    size = 16,
    face = rev(c('bold', 'plain', 'plain', 'plain','plain',
                 'bold', 'plain', 'plain', 'plain',
                 'bold', 'plain', 'plain', 'plain',
                 'bold', 'plain', 'plain', 'plain','plain', 'plain', 'plain',
                 'bold', 'plain', 'plain'))),
    plot.title = element_text(size=18, hjust = 0.5),
    legend.position = "bottom",
    legend.title = element_text(size=16), #change legend title font size
    legend.text = element_text(size=16),
    axis.title.x= element_text(size=16))+
  xlab("Marginal Means of Support for Candidate")



mmstable <- kable(mms)

save_kable(mmstable, "Tables/cand_mm.weighted.tex")


##MM heterogeneity
data_cand1 <- data_cand %>% filter(highsupportintervention==1)
data_cand0 <- data_cand %>% filter(highsupportintervention==0)

mms1 <- cj(data_cand1,   cand_choice~energy_policy + social_policy + migration_policy + cand_gender + cand_party, id = ~ResponseId,
           estimate = "mm")

mms0 <- cj(data_cand0,   cand_choice~energy_policy + social_policy + migration_policy + cand_gender + cand_party, id = ~ResponseId,
           estimate = "mm")

mmstablehigh <- kable(mms1)

save_kable(mmstablehigh, "Tables/cand_mm_het_high.tex")

mmstablelow <- kable(mms0)

save_kable(mmstablelow, "Tables/cand_mm_het_low.tex")


plot.data.mm1 <- create.plot.data.mm(mms1)
plot.data.mm0 <- create.plot.data.mm(mms0)



pl.mms.het<- ggplot(plot.data.mm1)+
  geom_point(aes(x = estimate, y = factor,  color = paste("High", " n = ", nrow(data_cand1)/10, sep=""), shape =  paste("High", " n = ", nrow(data_cand1)/10, sep="")))+
  geom_segment(aes(x=lower, xend = upper,
                   y = factor, yend = factor, color = paste("High", " n = ",nrow(data_cand1)/10, sep="")))+
  geom_point(data = plot.data.mm0, aes(x = estimate, y = factor, color = paste("Low", " n = ", nrow(data_cand0)/10, sep=""), shape = paste("Low", " n = ", nrow(data_cand0)/10, sep="")), position = position_nudge(y = -0.2))+
  geom_segment(data = plot.data.mm0, aes( x=lower, xend = upper,
                                          y = factor, yend = factor, color = paste("Low", " n = ", nrow(data_cand0)/10, sep="")), position = position_nudge(y = -0.2))  +
  geom_vline(xintercept = 0.5, lty = "dotted")+
  theme_bw()+
  scale_y_discrete(NULL, drop=FALSE,
                   labels = rev(c("Energy Policy",
                                  "Continue at Current Speed",
                                  "Slow Down Energy Transition",
                                  "Subsidies for Companies",
                                  "Subsidies for Low Income Households",
                                  "Migration Policy",
                                  "Refugee Camps Abroad",
                                  "Integrate Refugees Better",
                                  "Keep Refugees Away",
                                  "Social Policy",
                                  "No Pension Reform",
                                  "Early Retirement OK",
                                  "Raise Pension Age",
                                  "Party",
                                  "CDU/CSU",
                                  "AfD",
                                  "Bündnis90/Die Grünen",
                                  "Die Linke",
                                  "FDP",
                                  "SPD",
                                  "Gender",
                                  "Female",
                                  "Male")))+ 
  ggtitle("Marginal Means Candidate Conjoint,
by Belief in State Intervention")+
  theme(axis.text.y = element_text(
    size = 16,
    face = rev(c('bold', 'plain', 'plain', 'plain','plain',
                 'bold', 'plain', 'plain', 'plain',
                 'bold', 'plain', 'plain', 'plain',
                 'bold', 'plain', 'plain', 'plain','plain', 'plain', 'plain',
                 'bold', 'plain', 'plain'))),
    plot.title = element_text(size=18, hjust = 0.5),
    legend.position = "bottom",
    legend.title = element_text(size=16), #change legend title font size
    legend.text = element_text(size=16),
    axis.title.x= element_text(size=16))+
  xlab("Marginal Means of Support for Candidate")+
  guides(color=guide_legend(title="Belief in State Intervention"), 
         shape = guide_legend(title="Belief in State Intervention"))+
  scale_color_grey(start = 0.1, end = 0.5)



##MM difference
data_cand$highsupportintervention <- as.factor(data_cand$highsupportintervention)

mm_diff <- cj(data_cand,   cand_choice~energy_policy + social_policy + migration_policy + cand_gender + cand_party, id = ~ResponseId,
              estimate = "mm_diff", by = ~highsupportintervention)


mmstablediff <- kable(mm_diff)

save_kable(mmstablediff, "Tables/cand_mm_diff.weighted.tex")



plot.data.mm_diff <- create.plot.data.mm(mm_diff)



pl.mm_diff <- ggplot(plot.data.mm_diff)+
  geom_point(aes(x = estimate, y = factor))+
  geom_segment(aes(x=lower, xend = upper,
                   y = factor, yend = factor))+
  geom_vline(xintercept = 0, lty = "dotted")+
  theme_bw()+
  scale_y_discrete(NULL, drop=FALSE,
                   labels = rev(c("Energy Policy",
                                  "Continue at Current Speed",
                                  "Slow Down Energy Transition",
                                  "Subsidies for Companies",
                                  "Subsidies for Low Income Households",
                                  "Migration Policy",
                                  "Refugee Camps Abroad",
                                  "Integrate Refugees Better",
                                  "Keep Refugees Away",
                                  "Social Policy",
                                  "No Pension Reform",
                                  "Early Retirement OK",
                                  "Raise Pension Age",
                                  "Party",
                                  "CDU/CSU",
                                  "AfD",
                                  "Bündnis90/Die Grünen",
                                  "Die Linke",
                                  "FDP",
                                  "SPD",
                                  "Gender",
                                  "Female",
                                  "Male")))+ 
  ggtitle(paste("Difference in Marginal Means Candidate Conjoint,
by Belief in State Intervention"))+
  theme(axis.text.y = element_text(
    size = 16,
    face = rev(c('bold', 'plain', 'plain', 'plain','plain',
                 'bold', 'plain', 'plain', 'plain',
                 'bold', 'plain', 'plain', 'plain',
                 'bold', 'plain', 'plain', 'plain','plain', 'plain', 'plain',
                 'bold', 'plain', 'plain'))),
    plot.title = element_text(size=18, hjust = 0.5),
    legend.position = "bottom",
    legend.title = element_text(size=16), #change legend title font size
    legend.text = element_text(size=16),
    axis.title.x= element_text(size=16))+
  xlab("Difference in Marginal Means of Support for Candidate")


###merged figures for paper appendix

pl.mms.het.2<-  pl.mms.het + theme(axis.text.y=element_blank(),
                                   axis.title.x= element_text(size=12),
                                   plot.title = element_text(size=14, hjust = 0.5))

pl.mm_diff.2 <- pl.mm_diff + theme(axis.text.y=element_blank(),
                                   axis.title.x= element_text(size=12),
                                   plot.title = element_text(size=14, hjust = 0.5))


pl.mm_diff.2 <- pl.mm_diff + theme(
  axis.title.x= element_text(size=12),
  plot.title = element_text(size=14, hjust = 0.5))

mmggplotcand.all <- ggpubr::ggarrange(pl.mm, pl.mms.het.2, pl.mm_diff.2, nrow=1, common.legend = TRUE, legend="bottom", widths = c(0.8, 0.5,0.5))

ggsave(filename = "Graphs/candidate_conjoint.mmall.weighted.png",
       plot = mmggplotcand.all,
       width = 22, height = 10)


##########full table for split sample###########


split.highreg <- runreg_candidate_cov(data_cand, "highsupportintervention")
split.lowreg <- runreg_candidate_cov(data_cand, "lowsupportintervention")

split.highreg_a <- runreg_candidate_cov(data_cand[data_cand$candidateattention==1,], "highsupportintervention")
split.lowreg_a <- runreg_candidate_cov(data_cand[data_cand$candidateattention==1,], "lowsupportintervention")


row2 <- t(c("Demographic Controls",  "Yes",  "Yes", "Yes", "Yes")) %>% as.data.frame()

base_table <- modelsummary(list(split.highreg[[1]], split.lowreg[[1]], split.highreg_a[[1]], split.lowreg_a[[1]]),
                           stars=T, coef_rename = c("energy_policySlow Down Energy Transition"= "Energy Transition Slow Down",
                                                    "energy_policySubsidies for Competitiveness"= "Energy Transition Company Subsidies",
                                                    "energy_policySubsidies for Low Income Households"= 	"Energy Transition Low Income H. Subsidies",
                                                    "social_policyEarly Retirement OK"= 	"Early Retirement OK",
                                                    "social_policyRaise Pension Age"= "Raise Pension Age",
                                                    "migration_policyIntegrate Refugees Better"= "Integrate Refugees Better",
                                                    "migration_policyKeep Refugees Away" = "Keep Refugees Away",
                                                    "cand_gendermale" = "Candidate Male",
                                                    "cand_partyAfD" = "AfD",
                                                    "cand_partyBündnis90/ Grünen" = "Bündnis90/ Grünen",
                                                    "cand_partyDie Linke"= "Die Linke",
                                                    "cand_partyFDP" = "FDP",
                                                    "cand_partySPD" = "SPD"
                           ),
                           coef_omit = "region|age|income|gender",
                           add_rows = row2,
                           gof_omit = "AIC|BIC|Log.Lik|RMSE",
                           title="Split Sample Regressions, Candidate Conjoint Experiment",
                           output="latex")



tableheader.ts1 <- c(c(" " = 1,  header.tr = 4 ))

# set vector names 
names(tableheader.ts1) <- c(" ", header.tr)

high_a <- "Passed Attention"
all <- "All"

tableheader.ts2 <- c(" " = 1, all=2, high_a = 2)
names(tableheader.ts2) <- c(" ", all, high_a)

high_b <- "High Belief St. Int."
low_b <- "Low Belief St. Int."

tableheader.ts3 <- c(" " = 1, high_b =1, low_b = 1, high_b=1, low_b=1)
names(tableheader.ts3) <- c(" ", high_b, low_b, high_b, low_b)


regtable.b.c <- base_table %>% kableExtra::add_header_above(header=tableheader.ts3)%>% kableExtra::add_header_above(header=tableheader.ts2) %>% kableExtra::add_header_above(header=tableheader.ts1) 


save_kable(regtable.b.c , "Tables/candconj_split.weighted.tex")
